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In this paper, we study the collapse of a superhorizon-sized void in the 
O ' early, radiation-dominated universe using an improved general relativistic 



code. We find that in general a relativistic or nonrelativistic void collapses 
via a shock at the speed of light. This is true even if the outward velocity of 
the void wall is enormous. As the wall thickness decreases, the shock strength 
increases and the collapse time decreases up to a limit set by the shock tube 
solution. In addition, as the outward velocity of the wall increases, the col- 



^ . lapse time increases somewhat. When the collapse occurs in much less than 

a Hubble time outside the void, gravitational forces contribute negligibly to 
the solution at the collapse time; non-gravitational forces (caused by pressure 
and velocity gradients) shape the solution almost entirely. This is true even 
if the wall velocity is large enough that fiuid in the peak shocks outward sub- 
stantially during this time. Gravitational forces are expected to be dominant 
only after a void has collapsed. 
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I. Introduction 

Because the total energy density is of order tlie critical energy density today, if the 
early universe is described entirely by the big bang model, a severe fine-tuning problem 
exists. This problem is solved if inflation occurred in the early universe. First-order in- 
flation is one of the most interesting classes of inflationary models|]I|. Here, inflation ends 
in a particular region when a true- vacuum bubble is nucleated. Inflation ends everywhere 
when the universe is fllled with true vacuum bubbles. At this point, scalar fleld dynamics 
occurs which creates relativistic, fundamental particles. These particles thermalize much 
of the universe within a few Hubble times. Due to causality however, bubbles nucleated 
early on during inflation cannot be thermalized as quickly. Thus, since the inside of these 
bubbles are virtually empty, superhorizon-sized voids are created at the end of reheating.^] 
These voids are larger than the outside Hubble radius, cH~^^{ti}, (i.e. the Hubble radius 
outside the void), and can have enormous radii ~ 10'^'^ cH~^^{t{). The beginning of the 
radiation-dominated epoch then, is characterized by an approximately flat Friedmann- 
Robertson Walker (FRW) homogeneous and isotropic universe, punctuated infrequently 
by nearly empty, superhorizon-sized voids. 

The thermalization of the very large voids has been thought to take an enormous 
amount of time and therefore be in conflict with the tiny temperature fluctuations mea- 
sured in the microwave background 0. In a previous paper however (hereafter referred 
to as SV), it was found that these former estimates could be incorrect because the flrst- 
crossing time (the time taken for a photon in the void wall to reach the origin and thus 
make "contact" with a photon from the other side) was calculated incorrectly. In 0, it 
was found that the flrst-crossing time, At^, (i.e. the minimum thermalization time) is 
Ate/ H~^^^{ti) — (c~^-Rw/-f^(^t(^i))^- However, it was found in SV that the flrst-crossing 
time is actually much smaller than this: 



AtjH;^,{tO = <^,^'\c-'R^/H~^,iU)) [I + .5<l>,^-\c-'R^/H-^,{U))\ , (1.1) 

where the void's relative size is c~^i?w/-f^out(^i) < a^^^^, the log of its relative fluid po- 
tential is $in — a~^^^ 3> 1, a — Pin(^i)/Pout(^i) <^ 1, p is the energy density, and the 
subscripts "in" and "out" refer to quantities inside and outside the void, respectively. 
Because of the time dilation effect, the flrst-crossing time is smaller than the outside 
Hubble time if the void is empty enough. This effect occurs because the (relativistic or 
nonrelativistic) void has a large, negative fluid potential with respect to the relativis- 
tic background spacetime, which is caused by relativistic fluid effects. (In a matter- 
dominated universe, these fluid effects are tiny when the pressure is negligible.) Because 
the void is underdense, the time dilation effect is opposite of that for an overdense re- 
gion (e.g. a black hole) — time is dilated in an underdense region and contracted in an 
overdense region. If a superhorizon-sized void is empty enough then, it can collapse in 
less than an outside Hubble time. 

In this paper, we study in more detail the collapse of a superhorizon-sized void em- 
bedded in a FRW radiation-dominated universe. By collapse, we mean the evolution up 

''Following past convention, we loosely equate the Hubble radius with the horizon in the phrase 
"superhorizon-sized" . Horizon in this context is not to be confused with the particle horizon. 
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to the time the shock (created by the void wall caving in) collides with itself at the origin. 
In particular, we examine how large outward wall velocities affect the solution and also 
how gravitational forces affect the solution at the collapse time.0 We find that having a 
large outward wall velocity does not prevent a superhorizon-sized void from collapsing 
in less than an outside Hubble time. In addition, we show that if a void collapses in 
less than an outside Hubble time, gravitational effects are not important at the collapse 
time. This motivates us to compare the collapse of a superhorizon-sized void with that 
of a special relativistic void. We therefore derive the analytic collapse solutions in the 
slab limit, and test our code against these solutions. At the collapse time, thermalization 
and homogenization has not yet occurred, however. A companion paper will examine 
the subsequent evolution of a superhorizon-sized void after it has collapsed. 

The organization of this paper is as follows. In Section II we present the metric and 
equations of motion used to solve this problem. In Section HI, we describe the improved 
numerical scheme used to accurately find the solution near the origin, and we test this 
new code on the FRW radiation-dominated universe. In Section IV, we study the collapse 
of a superhorizon-sized void both analytically and numerically. In Section V we derive 
the solution for the collapse of an uncompensated void in the slab limit. Finally, Section 
VI contains a discussion of our results. 

II. Spherically Symmetric General Relativistic Fluid Equations 

We are interested in evolving a void embedded in a FRW universe. Because we know 
the initial conditions on an initial comoving time slice, we choose to work in Lagrangian 
synchronous coordinates. The metric then is 

ds"^ = -c2$2(t, r)dt^ + A'^{t, r)dr^ + R^{t, r)dn^, (2.1) 

where t is the coordinate time, r is the comoving radius, 2itR is the spacelike circum- 
ference of a sphere centered on the origin, and ln$ is the pressure gradient-induced 
"potential". We consider a perfect fluid with artificial viscosity. This fluid consists 
of particles with mass fi and temperature T. The stress-energy tensor in this case is 
rpap ^ c-2(p+p+Q)m°u'^ + (p+Q)^°^, where = (-c<l>-\ 0, 0, 0) is the fluid 4-velocity, 
p = nc^{l + e/c^) is the energy density, p = riT/p = ojne is the pressure, n is the mass 
density, e is the specific energy, Q is the artificial viscosity and w is a constant.0 The 
relativistic limit is obtained when e/c^ = T / (c^ufi) ^ 1. The viscosity is non-zero only in 
shocks. We assume that the total number of particles per comoving volume is constant, 
V ^{nu^) = 0, so that 

AtthR^R'/T = f{r), (2.2) 

where ' = d/dr, f{r) is chosen to be r^, |^ and 

F = R'/A. (2.3) 

■^Because we find in general that a void collapses at the speed of light, we will use the phrases "collapse 
time" and "first-crossing time" interchangeably. 
''Boltzmann's constant has been set to one. 



'^This choice only makes sense if F = 1 far outside the void, as will be clear following Eq. (2.15) 
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The fully general relativistic equations can then be written as0,0] 

R = (2.4) 

u _ ^(GnM AnGrAp + m\ T^Hp + Q)' 

^ ~ ^[ R^ ^ c2 J WnR' ^^-^^ 

M = -4:'k{p + Q)R^^U/c^ (2.6) 

^{p + Q){R^Uy 



_ _nHRHjy 

~ Rm' ^ ^ 



nR^R' 



(2.8) 



= -^^l±9l (2 9) 
where ' = d/dt, 

M' = Anc-^pR^R' (2.10) 

= l + iU/cf ~2GnM/{Rc^), (2.11) 

and W = 1 + [e + (p + Q) /n] /c^ (or c^nW = p + p + Q). Here, U is the fluid "velocity" , 
M is the "mass-energy" and W is the relativistic enthalpy. The energy density equation 
can be found by combining Eqs.(|2.7|) and 



p = -^p + p + Q)^-^. (2.12) 



Eq. (p.ll|) is the conservations of "energy" equation, with total, kinetic and potential 



"energies" of c^{T^ - l)/2, and -G^M/R, respec tively. For = 0, if a particle 

has velocity f in a Eulerian inertial frame, then T = 1/^1 — {vjcY and U = Vv (see Eq. 
( |A.1CI| )); r and U represent the two non-trivial components of the 4- velocity of the fluid. 
We can also write down an important auxiliary equation: 

f = -VU^ , (P + Q)' (2.13) 
{p + P + Q)R' ^ ' 

Thus the "energy" c^(r^ — l)/2 can only decrease in the presence of pressure gradients. 
In SV, Eq. ( |2.13| ) was used to argue that if F 1 in the void wall, the wall might slow 
down in a relatively short amount of time. We will find later that this is true in some 
cases. 

It is necessary to have artificial viscosity in the code due to the presence of shocks — it 
is added to physically alter the otherwise incorrect solution. Artificial viscosity dissipates 
just enough energy over a shock to satisfy the (exact) jump conditions, which cannot occur 
in a perfect fluid. These conditions are derived in the appendix, and a demonstration of 
the codes ability to give the correct solution with the following form for Q will be given 
there (Figure 13) and later on in this paper (Figure 9). 
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The form of the artificial viscosity used here is a generahzation of the nonrelativistic 
form0, and allows relativistic shocks to penetrate into relativistic fluids. It does not, 
however, alter the solution substantially in non-shock areas. 

Q = k'^nWT-^iU'fdr'^ for U' < 

Q = otherwise. (2.14) 

The constant k"^ is of order one, and is proportional to the number of grid points in the 
shock. For relativistic shocks penetrating non-relativistic fluids, Eq. ( p.l4| ) becomes the 



expression obtained previously . The expression from this reference however, does not 
give enough viscosity over strong shocks when the fluid in front of the shock is relativistic, 
because Q <^ p. (It is necessary for the artificial viscosity to be of order the pressure 
in the shock region, in order to dissipate enough energy). Eq. ( p.l4| ) does, however. 
Assuming that the fluid in front of the shock is stationary, then k'^AU'^ ~ — 1 over 
the shock front, where we have used the shock jump conditions Eq s. (|A . 2 6|) - (|A . 2 7|) in the 



strong shock limit. Then Q ~ {p + p)U^ /(T^c^) ~ p(F^ — 1)/F^ ^ p ~ p, the desired 
result. 

If the fluid is homogeneous and isotropic, Eq. ( |2.2| )- (|2TT|) reduce to the FRW equa- 
tions, with 



R = ra, $=1, T = ^1- Kr^/c'^ = ^l + {RH/c)'^{l-n), 
M = Anc-^pRyS, H = $"^a/a = U/R, (2.15) 

where a{t) is the cosmic scale factor, H is Bubble's "constant", f2 = 1 + n/{HaY = 
SnGNp/i'ic^H^), f = r^/^l — nr'^/c^, and k is —1, or 1 for negative, zero and positive 
spatial curvature, respectively. In a spatially flat (k = 0) FRW universe then, U = 



c~^RJ'&ttGnp/^ = / R. Note that F(t,r) = 1 in a spatially flat universe, even 

though the fluid everywhere is moving. 

III. Numerical Scheme and Tests of Improved Code 

In this section we describe the improved numerical method used to solve the equa- 
tions given in the previous section. We also test the code on the FRW radiation- 
dominated solution. Readers more interested in the physical results should skip to section 
inc. The code is greatly improved over that used in SV, because we have now imple- 
mented numerical regularization into the difference equations This is done to get 
rid of instabilities that develop at the origin from the inaccurate calculation of e, h and 
f^FLUiD = — F^$(p + QY / {WnR') there. As we will see in a moment, the difference equa- 
tions used in SV were only accurate to (AR/R)'^, where AR is the grid size, which is 
or order one at the first few grid points surrounding the origin. Thus, no matter how 
small AR is chosen to be, e, n and f/pLuiD could not be calculated correctly for these grid 
points. It is true that decreasing AR confines the inaccuracies to a much smaller physical 
volume. The problem, however, is that the inaccuracies can become instable when violent 
fluid processes take place at or near the origin. An example of such a process is when 

■^It turns out that scheme 1 from SV for the h equation is differenced correctly, however. 
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a spherical shock colhdes at the origin. These instabihties then cause the code to fail. 
Numerical regularization ensures that the error made in calculating e, n and t/pLuiD is of 
order (Ai?//)^, where I is the characteristic length scale of the problem. With numerical 
regularization then, the error made can be a small as you like by choosing AR sufficiently 
small. 

A. Numerical Regularization 

We start by differencing the derivative term in Eq. (^) and (^Tsf) , (R'^U)' /{R'^R'), 
at spatial grid point j as A{R^U)/{R^AR), as was done for the n equation in scheme 2 
of SV. We approximate R and U at adjacent grid points j ± 1 by Taylor expanding: 



Rj±i = Rj ± Ar{drR)j 
Uj±i = Uj±Ar{drU)j. 

If we forward difference this derivative term at grid point j, we find 



(3.1) 



-Rj_i_i(-Rj+i — Rj) Rj{drR)j 



+ 



R]m), 



2UjRj{drR)j |l + 



Arid^] 
2R, I 

{drR)] 



(3.2) 



Backward differencing is obtained by substituting Ar — —Ar into Eq. ( |3.2| ). The forward 
(backward) differencing gives the predicted (corrected) solution, as will be described in 
section B. The predicted and corrected solutions are then averaged together. When this 
is done, the error made in calculating e or n at the new time step is approximately 
{AR/ RjY . This truncation error is small for all but the first few grid points near the 
origin. There, the truncation error is of order one, so that inaccurate calculation of e 
and n will result when U' ^ 0.| If our grid is equally spaced in R so that Ri = AR/2, 
i?2 = 3Ai?/2, etc, then for j = 1, 2, 3 and 4, the truncation error in the R'^{drU)j term is 
about 400%, 44% 16% and 8%, respectively. A similar truncation error problem occurs 
when differencing f/pLuiD as in SV: f/pLuiD oc A(p + Q)/A{r^). 

We can rectify these problems using numerical regularization. Near the origin, we 
expand p, Q, R and U in a Taylor series in r. Then we use symmetry arguments to deduce 
the actual form of each solution there. Using the fact that p{—r) = p{r), Q{—r) = Q{r), 
R(-r) = -R{r) and U{-r) = -f/(r) along with the fact that R{0) = f/(0) = 0, p(0) ^ 
and Qj oc {U'Y, the lowest-order terms in the series at grid point j are 



Pj 
Rj 



Po+P2r^ + Oirj) 
Rirj + R^r] + 0{r]) 
Uir, + Usr^ + 0(r|) 

= Q2 [f/i(r,+i - n) + f/3(ri+i - rl) + 0{rl rl^,] 



Q2iu, 



k+l 



, (3.3) 



^In SV, e was determined via calculating A(i?^C/)/(r^Ar). It can be easily shown that this gives a 
truncation error of order 0{{AR/ Rj)"^) also. 
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where the coefficients po, p2, etc, are constants, and where k = j and k = j — 1 for the 
forward and backward derivatives, respectively. Because (i?^f/)j+i = RlUir^_^_^ to lowest 
order, it is clear that we should instead difference (R'^U)' /{R'^R') as 



{R^uy 
Rm' 



R3 



Rj 



(3.4) 



for the forward derivative, and similarly for the backward derivative. We can verify this 
by plugging Eq. (|]3D into Eq. ( ^ : 



,-R|+if/j+i 



R3 



i?3 



3Ui 



1 + (f/3/f/i + 2Rs/R^)r' 




Ar / Tj 



1 



(1 + Ar/r-j f - 1 



1 + {m^/Ri)r 




Ar/vj)^ — 1 
Ar/rjY — 1 



-I -1 



(3.5) 



The correct result from Eq. ( ^.3] ) is 3Ui/Ri. Since the quantity in curly brackets is 
of order 4 or less, U3/U1 ~ 0{{drRp-'^) and i?3/i?i ~ 0{{drRp~^), where / is the 
characteristic length scale of the problem, the truncation error is of order [AR/lY for 
Rj ~ AR. Therefore, the error made in calculating this derivative can be made as small 
as desired by decreasing the grid size AR. 

Now we turn to the derivative term in (J-piAjm'- {p + Q)' / R' ■ Using Eq. (|3.3| ), we 
see that the obvious way to difference p' / R! is 2Rj{pj+i — pj)/[{Rj^i + Rj){Rj+i — Rj)] 
for the forward derivative, and similarly for the backward one. However, this causes 
a problem at j = 1 for the backward derivative because po = pi and Rq = —Ri, but 
Ri{pi ~ Po)/[{Ri + Ro){Ri — Ro)] 7^ 0. Instead, we use the simple difference scheme 



P. 
R' 



Pj+i - Pi 
Rj+i - Rj 



(3.6) 



for the forward derivative, and similarly for the backward one. Plugging Eq. (|3.3| ) into 
the previous equation, we get {2p2rj/Ri){l + Ar/{2rj)) and {2p2rj/R\){l — Ar/{2rj)) 
for the forward and backward derivatives, respectively. Upon averaging the forward and 
backward derivatives, we get the correct result, 2p2rj/Ri, to lowest order Note that in 
this case, the backward derivative at j = 1 is zero, thereby avoiding the problem caused 
by the former difference scheme. 

We also difference the Q' / R! term in the same way as Eq. (|3.(j|) with p replaced by Q. 



Using Qj from Eq. 

Qj+i ~ Qj 



with k = j, the forward derivative is 



l2Q2UiU3Ar\ 



Rj+i - Rj 



Ri 



2-3 



+ 



36Q2UiAr'r 



3Ar 10 ( Ar 




(3.7) 



''Although pi and Ri will be slightly different for the forward and backward derivatives, their varia- 
tions are of higher order. 
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The backward derivative can be obtained by substituting — Ar in for Ar in the previous 
expression. We can now average the forward and backward expressions and obtain to 
lowest order0 

( in / 2' 

1 + (3f/3/f/i)r2 



12Q2UiUsAr^rj 



Ri 




(3.8) 



The correct result, from Eq. (^, is 12Q2f/3ArVj(f/i + 3f/3r|)//?i. But because U3/U1 ~ 
0{{drR)^l^^), where / is the characteristic length scale of the problem, the quantity in 
brackets is [1 + 0{AR/l)^] for Rj ~ AR. Therefore the error can be made as small as 
desired by decreasing the grid size AR. 

B. Numerical Procedures 
We set up a grid with Jb + 1 grid points which is equally spaced in R such that 
Ai?(ti) = i?(ti)j+i - R{U)j for j e [0,jb], with Rq = -AR{U)/2 and Ri = AR{U)/2. 
Then we specify p(i?(ti)) (or M{R{ti))) and T{R{U)) (or U{R{t-^)). When T{tuR) is 
given, we always choose to have the fluid moving outward initially, ?7 > 0, using Eq. 
( p. lip . Assume that the fluid is initially an isentrope, so that the specific energy, e, can 
be determined via solving 

fl) = , '^]-"^\,, { (1 + ^(*.. Hyc') (3.9) 



l+e((|,fifl)/c2 \4k,RB) 

iteratively, where the subscript "5" denotes the quantity at the outer boundary. In 
addition, $(ti, Rb) = 1, and the artificial viscosity is initially zero: Q{ti, R) = 0. Finally, 
Tj and $(ti,i?) can be determined by integrating Eqs. (|2.2|) and ( p.9|) from R = and 
R = Rb, respectively. It is important to note that for p = p/3 and (5 = 0,$ can be 
solved exactly on the initial time slice: 

$(ti,i?) = <l>out(ti) [j^) ' e(ti,i?)/c2 > 1, (3.10) 

where "out" refers to quantities outside the void at j = jb- We define p as the energy 
density where the fluid becomes nonrelativistic: zjc? ~ 1. Since the contribution to $ is 
negligible when the fluid is nonrelativistic, 

$(ti, R) = $out(ti) (^^) ^ for p(ti, R) < p. (3.11) 

Inside the void, p = pin(ti) for a relativistic void, and p ~ g{c^pY fo^ ^ nonrelativistic 
void, where g is the number of degrees of freedom, p is the particle mass, and "in" refers 
to quantities inside the void. For a nonrelativistic void in a relativistic background fluid 
then, 

<fin(ti) ^<fout(ti)eout(ti)/c^ (3.12) 



* Again, we ignore the variations of ?7i, C/3, Q2 and i?i in the averaging, since they are of higher order. 
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We use the MacCormack predictor-corrector method to integrate the equations of 
motion. Suppose we know all quantities on the i^^ time slice for all spatial points j. 
Using Eq. ( |2.7| ) as an example, we first predict the new quantities (with forward dif- 



ferencing) for all j: 



p j 



At 3n]^] 



{R 



'j+i Uj+i 



Rr 



)■ 



After using similar methods to obtain f/p*^^, Rp^^^ ■, Mp^^^ and Cp*^^, we integrate Eq. 



^ inwards from the outer boundary using the 4*''-order Runge-Kutta method with 
linear interpolations to determine ^p*"^^- The "predicted" viscosity (to be used in the 
corrector step) is: Qp'^^ = k'^niW^iTi)-^{U}^^ - U^f from Eq. (CT ). We then in- 
tegrate again (with backward differencing) and average to obtain the corrected values 



n 



= .5(n} + n;f - At Sn^^^pJ (R;/ Up/ - Rp)__, Up;__^/{Rpf - Rp;_, The 



other quantities are obtained similarly. (Note that the velocity is backward differenced 
also when calculating Q^^^)- After the n*^ corrector step, we calculate the time step 
for the {n + l)*'^ integration, (At)"^^, which must satisfy the Courant condition and the 
condition that n, e and M change slowly enough at each grid point. (This last condition 
is necessary for gravitational expansion and contraction). Thus 



(At) 



n+l 



mm 



/ 



V 



TDn 



cs) 



f f — f — - 

Iri" I ' I ^'^ I ' 



(3.13) 



where {cs)) = y + ^){p] + Qj)/{^j^j) is the speed of sound,[] / < 1 is a constant, 
and n/\n\, etc, are averaged over a few surrounding grid points. 

At the origin, reflecting boundary conditions are used. And at the outer boundary. 



^{t,R 



B) 



1, p' 



n 



and Q{t,RB) = 0, so that Eqs. 



can be 



easily integrated using the MacCormack method. Then n and e are determined using the 



'free-string" condition: rij 



3b 



n 



'IB- 



and e 



Jb 



C. Void Initial Conditions 

A void is defined by its radius wall thickness AR^, and its relative energy density 
a ~ Pin(ti)/ Pout(ii) < 1- If c~^Rw/H~^^{ti) > 1, the void is said to be superhorizon-sized. 
The void wall can be either compensated or uncompensated. If the void is initially 
uncompensated, then the energy density distribution used here is 



p(ti, R) = .5pout(^i)[(l + tanhx) + a{l — tanhx)]. 



(3.14) 



where x = {R — R^)/ AR„. Thus the energy density increases monotonically to its 
outside FRW value, Pout(''^i)- If the void is initially compensated, then the "mass-energy" 
distribution used here is 



M{ti,R) = .5c"^47rpout(ti) [(1 + tanhx) + a(l - tanhx)] /?^(ti)/3. 



(3.15) 



■'We include Q in the expression for {csYj because we want the effective pressure. In practice, it 
makes httle difference. 
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In this case, M reaches its FRW value outside the void. Therefore, the excess mass- 
energy density in the void wall compensates for that missing from the void. Using Eq. 
( p.lOp , one finds that the energy density for a compensated void is 

P = Pout(ii) \ ^ . „ ^ ,2 + i [(1 + tanhx) + a{l - tanhx)] I . (3.16) 
[ 6Ait„ cosh X 2 J 

The actual wall thickness can be found by determining the radius, -Rinner, at which the 
energy density is approximately that inside the void. Setting p(ti, i?inner) — 3pin(ti) and 
using Eq. ( [j.l5| ), we can solve for -Rinner : 



-Rinner -Rw , /3A-Rw -i r> i ,1 /3A-R^ , , 

AP -^"^V^ a ~ 1.21ogio« + lnW— , (3.17) 

^-iCw V -'dinner V -'''inner 

where we have assumed that 3A-R„ -C -Rw Thus the actual wall thickness is larger than 
A-Rw, since it is of order or greater than | log^g a| A-R^. 
The initial velocity distribution is determined via 

r(ti, i?) = 1 + (r„ - 1) exp[-(i? - R^ f/i24)], (3.18) 

where cxr and > are constants. For this function, F is equal to one inside and 
outside the void. In the void wall, F can be greater than or equal to one. 

In many instances, we would like to match the position of the inner edge of the 
energy density distribution with that of the inner edge of the velocity distribution. At 
-R = -Rinner, the excess kinetic energy should be no larger than the "potential energy". 
Setting F = 1 + /5 with /3 < 1 and using Eq. ( pTl| ), we require 20^/3 < 2GnM/R. This 
occurs when 

where b < 1. Using Eq. (|3.18| ), the width of the velocity distribution then is 



ar = , . (3.20) 

^21n[(F„-l)//5] 

For nearly all of the simulations in this paper, we set 6 = 1 to determine cr- As an 
example, if -Rw = 50, A-R^ = 1, a = 10~^, F„ = 6, H~^^{ti) = 2 and b = 1, then 
i?inner = 44, /? ~ .024, and ar ^ 1.8. 

As a final note, a void from first-order inflation is compensated and has a large 
outward wall momentum|P]. The mass (M) outside the void then, is the same as if the 
void were not present, and the wall "velocity" is large: F„ > 1. 

D. FRW radiation-dominated model 

In this section, we test our improved code out on the FRW k, = radiation-dominated 
model. We also apply a convergence test since the exact solution is known. The relative 
error in q, where q denotes any quantity, is defined to be Cj = (gj — q{ri))/q{ri), where 
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g(rj) is the exact solution and qi is the numerical solution. We obtain a global measure 
of the error by defining 



1 ^ 



Ll = J^ll\^^l (3.21) 



=1 



where is the total number of grid points. This error is proportional to the grid spacing 
to some power: Li oc Ai?'^, where s is the convergence rate. If s ^ 2, the code is 
second-order, as desired. These tools have been used previously to test codes in other 
applications |1T0[|. 



For the simulations in this section, we set a = 1, Gn = 1, c = 1, o; = 1/3, t\ = 1, 
eout(^i)/c^ = 10^ and k"^ = 0. We examine a fiuid which is initially relativistic, ho- 
mogeneous and isotropic, with T{ti,R) = 1. The solution is R{t,r) = R{ti,r){t/tij^/'^ , 
47rGjvPhom(t) = 3cV(8t2) and U = R/{2t). 

We set c-^RB/H-^{ti) = 250, Ai?(ti) = 2.5 and C = .3, and run the code to t = 1.1. 
These are the same initial conditions used to make Figure 4 in SV. Figure 1 shows the 
relative error in p plotted versus the scaled radius at this time for / = .01, .005 and .0025 
as solid, dotted and dashed lines, respectively. We see that the relative error in p is not 
only independent of i?, but is also very small, as compared to Figure 4 in SV. For instance, 
for / = .01, the inner and outer boundaries have —6% and —2% errors, respectively, while 
it is only 3 x 10~^% here. Thus, the energy density is no longer underpredicted for the 
first few grid points near the origin, which is an important consequence of numerical 
regular ization. 

We now calculate the convergence rate in this model. Since Gjy 0, the size 
of the time step is determined by calculating the smallest of the following: Atcg = 

CAR/iT^cs) - V^CAR{ti)^/t/tiand At„ = J \n/n\ ~ / 2t/(3ti) from Eq. ( pTT^ ). Thus 
we will test our code in two limits: when Atc^ < Atn (so that decreasing AR(t{) decreases 
At), and when Atn < Atcg (so that decreasing / decreases At). Running our code in the 
former limit, we set Rb = 40, / = .5 and C = .001, and run the code tot = 1.035. (Note 
that we have made / large and C very small in order to achieve this limit). In Table 1, 
we show the values of Li{q) for q = n, R, U, M and e for various AR(ti), and calculate 
the convergence rate. We see that the rate is roughly 2 for all variables, although there 
is some variation. 

We also run our code in the latter limit, and thus set Rb = 40, AR(ti) = 5.0 and 
C = .3, and run the code tot = 1.014. In Table 2, we show the values of Li{q) for various 
/. We also calculate the rate s, defined by Li oc f^, in order to determine the speed with 
which the code converges when the time step is determined solely by gravity: At„ oc /. 
Again, the rate is roughly 2. 

IV. Collapse of a Superhorizon-sized Void 

In SV, we found that the superhorizon-sized voids examined collapsed at the speed 
of light. The collapse occurs because the fiuid in the void wall shocks inward from the 
wall pressure and velocity gradients, and it does so at the speed of light because the 
shock is in general strong (see Eqs. ([A.24| ) and ([A.30| )). In addition, it was found that if 



the relative energy density in the void is small enough, the collapse occurs in less than 
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an outside Hubble time (i.e. a Hubble time outside the void). We ran many examples 
for which the wall energy, (F^ — l)c^/2, was zero. However, we only ran one example 
{Figure 11) for which it was non-zero and large. The void still collapsed in this case. 
At the collapse time, Atcoiiapse/-f^(^t('^i) — ■'^^ there was substantial pressure and velocity 
distortion at the peak and in the peak area. We would like to understand this example 
and generalize to other situations. Do all superhorizon-sized voids collapse, even if the 
wall energy is enormous? What causes the distortion in the peak area, and does it always 
occur? And how does the wall energy affect the collapse time? 

In order to answer these questions, we will take two approaches. First, we calculate 
analytically the time scales for change of the fluid in the void wall as a function of the void 
radius, wall thickness, wall velocity and relative energy density, a. Second, we answer 
these questions by showing the results of many numerical simulations. We consider 
different evolutionary outcomes caused by changing the void's inside energy density, wall 
thickness and excess wall energy. We then determine which forces (gravitational or non- 
gravitational) are responsible for accelerating the fluid into the void and for distorting the 
initial peak and peak area distributions by comparing runs of superhorizon-sized voids 
with those where gravitational forces are neglected. It turns out that gravitational forces 
are negligible if the collapse occurs in much less than an outside Hubble time. This may 
seem confusing, because superhorizon-sized voids obviously evolve under gravitational 
forces. However, the time scales for which gravitational forces substantially change the 
solution can be much larger than that for which fluid forces (gradient forces) change the 
solution. When this is the case, gravitational forces negligibly affect the solution for some 
time. Then, the conclusion that non-gravitational forces have caused virtually all of the 
observed changes in the void is a correct and important result. 

A. Time Scales for Fluid Motion in the Void Wall 

In section B, we will calculate numerically the evolution of many different superhorizon- 
sized voids with different initial conditions. From these simulations, detailed evolutionary 
results can be obtained. However, it is difficult to generahze these results to situations 
with different initial conditions. For this reason, we calculate analytically the time scales 
on which the fluid in different regions of the void wall changes. This will allow us to gain 
intuition over the entire parameter regime of initial conditions. 

Consider a compensated superhorizon-sized void with initial energy density given by 
Eq. (|3.16| ) with a ^ 1 and 3AR^ <^ R^.^ In addition, the initial velocity is determined 



via Eq. (|3.18|) , and the pressure is p = up, where a; = 1/3 (a; = 0) if the fluid is relativistic 
(nonrelativistic). We assume that the fluid outside the void and in the peak area of the 
wall is relativistic, whereas the fluid in the void may be relativistic or nonrelativistic. 
The wall pressure, "mass" and "velocity" at the peak, R = R^, is then 

Pw((,) ^ (4.1) 
MAt.) - \('-^rA (4.2) 



fe 



This is not a requirement that the void wall be very thin, because the actual wall thickness is greater 



than or of order | logj^Q a| Ai?„ from Eq. (3.17) 
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Note that the wall energy, c'^(T^'^ — l)/2, is the major contributor to the wall velocity 
only if r^>c~'R^/H-;,{U). 

We now calculate the dynamical time scale on which the fluid at a particular position 
changes. To do this, we first calculate the time for each of the fluid variables R, U, M, p 
and r to change by assuming that all variables are constant and are given by their initial 
values.! The change in velocity is the sum of the gravitational (Gat 7^ and (p + Q)' = 0) 
and "fluid" {Gn = and {p + QY 7^ 0) contributions, which we consider separately. The 
fluid variable which changes the fastest then sets the time scale for change of all variables. 

As an example, we calculate the time scales for which p and t/pLuiD change. We 
set the quantities in Eqs.( |2.5|) and ( |2.12| ) equal to their initial values, and calculate the 



times at which At\p\/p ~ 1 and Atlt/pLuiDl/f^ — 1- The time scales for change of p and 
U FLU ID then are 

respectively. A similar procedure can be followed to obtain At[R], At[UGRAv], At[M] 
and At[r] from Eqs.( p.4|) , ( p.5|) , (2.6) and (|2.13|) , respectively. The time scale for change 



of p, for example, is given by the smallest At obtained, rather than At [p\. This is 
because we can no longer assume that the fluid quantities are equal to their initial values 
if some other quantity is changing on faster time scales. If U flu id changes the fastest, 
for example, R and U (and therefore p) will also change on this shorter time scale. 

We are interested in the following two regions: the fluid near the peak (i? = i?„±Ai?w) 
and at the peak {R = i?„) of the wall. When x = —1, the energy density is p = .42pw ^ 
Pout(^i) and the spatial derivatives are (R^U)' /{R^R') ~ U^/AR^ and p'/R' ~ p^/AR^. 
If the void is larger than four times the outside Hubble radius, then the energy density 
changes the fastest 

1/4 



1 , Ai?„ [R^ARj 

Atpeafc area = At[p] = ^—$-1 = .39 ^ 



/AD X 3/4 

< -ssi/o-Uti) ( =F ) > (4-5) 



r\-/ 



V Rv. 



where we have used the fact that > R^/ {^/2H~^^(ti)). The same time scale for change 
is obtained at a; = 1 as long as Ux=i <^ f/w. As the wall gets thinner or the wall energy 
gets larger, the time scale for change decreases because of the steeper velocity gradient. 

Thus we find the interesting result that the fluid in the peak area always moves in 
less than an outside Hubble time. Note that the time scale for change is independent of 
a ~ Pin(^i)/Pout(^i)5 which is not very surprising because the velocity gradient will hardly 
increase if Pin(ti) is decreased. 

'The time scales of change for n and e differ from that for p by a number of order one. 
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The second region of interest is at the peak {x = 0), where p' = T' = 0. Then 
{R^Uy /{Bj^R!) ~ 2f/„/i?^ + (i?^iJout(ii))V(4f^wAi?„), so that mass-energy sets the time 
scale for change: 

Atpeak = At[M] = 3$w -TT- - 1-9^ 7T 

/AD X 3/4 

< (^) , (4.6) 

which is 5 times larger than Atpeak area- Again the time scale for change is less than the 
outside Hubble time. Therefore a superhorizon-sized void must collapse in at least less 
than a Hubble time if the fluid in the peak area is to remain virtually stationary during 
the collapse. 

Everywhere in the void wall, the time scale for change of the energy density is At[p] ~ 
(1 + ct;)^^$^^A_R„/[/„, which is greater than or equal to the time scale of change at that 
position. As we move from the peak area into the void, decreases because the fluid 
"potential" decreases due to the time dilation effect. Thus the time scale for change 
decreases as R decreases, so that the fluid near the base of the void accelerates into the 
void on time scales faster than fluid in the peak area can move; a superhorizon-sized void 
cannot prevent its own collapse by having excess energy in the peak wall area. 

We can compare the time scale for change of a superhorizon-sized void to a void 
with nearly identical initial conditions, but for which the gravitational forces are zero:Q 
Gtv = 0. This is a SR void. At x = — 1, we assume that r„ > F > 1. Then At^pi^uij^ is 
the smallest only if F — 1 <C 1 and F — 1 < .1/ (F„^ ~ !)• Otherwise, At[p] is the smallest: 

At,,akarea = At[p] ^ -^$"1^ = .39^ (4.7) 

When Ux=i <^ the time scale for change is given by Eq. ( ^?7| ) also. This agrees 
with the results in the GR case when F^ ^ i?„/(cif(^t(ti)). When F„ < R„/{cH~^^{ti)), 
Atpeak area IS Smaller in the GR case because is larger. 

The time scale for change at the peak is quite different in the SR case however, 
because the energy density now sets the time scale for change: 

Q p p 5/4 

AW-A«W^3*;'--.24^^. (4.8) 

This is larger than the GR result by the factor i?„/Ai?„ ^ 1. Thus, it takes much longer 
for the fluid at the peak to move for a SR void as compared to a GR void. This has 
interesting consequences. Because the energy density of fluid at x = — 1 decreases more 
rapidly than that in the peak, the pressure gradient in the upper wall area increases. 
This causes F„ to decrease, which slows the wall down rather quickly. We will see this 
in section Bl numerically. 



'R^, Ai?w, ri{ti,R), e{ti,R), <i>{ti,R), R{ti,r) and T{ti,R) are the same. Only U{t;,R) differs 



through Eq. (2.11) 
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As an example, suppose we consider a super horizon- sized void with = 50, A_Rw = 
1, r„ = 6 and a = 10~^ embedded in a FRW k = radiation dominated universe with 
-^oiit(^i) — 2- Then U„ = 19 so that the time for the peak and peak area energy density to 
change is .13H~^^^{t{) = .27 and .028if^t(ti) = .055, respectively. If we now consider a SR 
void with the same parameters and for which F > 1.003, then f/^ = 5.9 so that the time 
scale for the peak and peak area energy density to change is 5.4 and .18, respectively, 
which are much larger times. 

It is interesting to calculate the initial conditions a superhorizon-sized void must have 
in order for its peak area energy density to remain roughly constant at the collapse time. 



From Eq. (4^), it is clear that the collapse must occur in less than an outside Hubble 



time. From Eq. ( |1.1D , this occurs if $in c ^R„/Ho^^{ti) < 1, where $in is assumed to 
be nearly constant. In this case, the collapse time is 



At, 



collapse 



(4.9) 



In general, $in is larger than its initial, isentropic value, because the fluid is not isentropic 
across the inbound shock which develops. Assume that $in = (Pout(^i)/p)^^^, where p is 



defined after Eq. ( |3.11| ) and p > l.Q Letting AtcoUapsc < Atpeak area, we find that 



Poutiti] 



3/p 




-2/p 



(4.10) 



Therefore, if p/pout(^i) is small enough, then at the first- crossing time the fluid in the 
wall will not have moved very much. It is important to emphasize that if AR„/R„ -C 1 
or f/„ ~ cF^ 3> 1, then the collapse must occur in much less than an outside Hubble 
time in order for this to occur. 

B. Collapse — Numerical Results 
In sections Bl and B3, we compare the collapse of a superhorizon-sized void with 
that of a void with similar initial conditions, but for which we neglect gravity. (This 
latter void is a special relativistic (SR) void). The purpose of this is to compare which 
evolutionary changes in the void are due to gravitational and non-gravitational (i.e. fluid) 
forces, since any effects which appear similarly for both SR and GR voids cannot be due 
to gravitational forces. We find that if a superhorizon-sized void collapses in less than the 
time scale for change in the peak area, the solution can be approximated very well by the 
solution to a SR void with the same initial conditions. For these voids then, gravitational 
effects are not important during the collapse. In addition, we find that even if the collapse 
time is larger than the time scale for change in the peak long as the collapse occurs 

in much less than an outside Hubble time, gravitational effects are unimportant. These 
numerical results are important, since they show that the gravitational contribution to 
the solution is negligible at the collapse time if a superhorizon-sized void collapses in 
much less than an outside Hubble time, regardless of the details in the wall area (e.g. 
wall thickness, wall energy, etc). 



"From Eq. (3.11), p — I only in the isentropic limit. As an example, p = .366/. 25 = 1.46 for much of 



the collapse of a thin- walled imcompensated relativistic void, using Eq. (5.24) 
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For Figures 2-7 and Tables 3-4, we choose T(ti,R) and p{ti,R) from Eqs.( |3.18| ) and 
( p.l6| ), respectively. We then integrate the equations of motion twice, starting with the 
same initial conditions.^ The first time we set Gat = 1 (the "GR" case), which gives the 
solution for the evolution of a super horizon- sized void. The second time we set Gat = 
(the "SR" case), which gives the solution for the evolution of a special relativistic void. 
The pressure in the GR (SR) case at time t is defined to be pcRit, Ri) {psR(t, Ri)), where 
Ri = R{ti,r). 

In all figures which follow in section IVB, we set C = .3, c = 1, u = 1/3, / = .01, 
ti = 1, eout(^i)/c^ = 10^ and 47rpout(ii) = 3/8. Thus the initial Hubble time outside the 
void is if^t(ti) = 2ti = 2, and all voids are relativistic.0 
1. Gravitational Effects during the Collapse of a Superhorizon-sized Void 

We start with compensated voids for which R^ = 50, AR^ = 1, T(ti,R) = 1, and 
AR{ti) = .25. Note that the GR void is superhorizon-sized, since R^ / H~^^{t\) = 25. 
For a = 10"^, 10"'^ and 10~^°, we set k"^ = 1.5, 1.7 and 1.7, respectively, and we run our 
code twice to t = 4.1, t = 1.5, and t = 1.08, respectively. 

In Figure 2, we show the fractional difference in the pressure, [pcRit, Ri) —PsR(t, Ri)]/ 
PGRit, Ri), as a function of the initial radius. For a = 10"*^, 10"'' and 10~^°, we plot solid 
triangles, open boxes, and a dashed line, respectively. As the inside energy density 
decreases, the difference between the SR and GR solutions decreases also. For a = 10~^, 
10~^ and 10~^°, the absolute value of this difference is at most 350%, 110% and 20%, 
respectively. A similar result is obtained for an initially uncompensated void. 

When a = 10~^, the solutions at t = 4.1 are very different. Most of the fluid not 
directly behind the shock is moving outward in the GR opposed to the SR 

case. In addition, the large void wall pressure, p„, has redshifted away in the GR case, 
causing the pressure gradient to be almost nonexistent in the former wall region. This 
is in contrast to the SR case, which has hardly changed at all in the peak area. The 
fact that the solution at t = 4.1 has changed considerably in the GR case is not very 
surprising, because the time scales for change in the peak and peak area are .29 and .059, 
respectively, from Eqs. ([4.6|) and ([4.5|) . 

In Figure 3, we show voids identical to those of Figure 2 but for which = 6 and 
(Tr = Ai?„. For this value of r„ however, U„/c = 19 > r„, so that the velocity is about 
5% higher than for Figure 2. Again, as the GR void becomes emptier, the solution looks 
more like the SR solution. For a = 10"^'', gravitational effects can be neglected during 
the collapse of the superhorizon-sized void. Note that the relative pressure difference 
for a = 10""^ is only 50% in the shock opposed to 300% in Figure 2. This is 

because the SR shock here is weaker since more wall fluid is pushed out with the larger 
wall velocity. A similar result is again obtained for the uncompensated case. 

In Figure 4, we plot the radii of seven comoving observers as a function of time for 
the voids in Figure 5" with a = 10~^ and a = 10~^°. The solid and dashed lines represent 
the GR and SR results, respectively. The figures show that when a = 10""^, the SR 

°The only variable which differs between the GR and SR cases initially is the velocity U{ti,R). All 
other variables are identical. 

PA relativistic (non-relativistic) void is when the fluid in the void is relativistic (non-relativistic). 
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and GR results deviate quickly. Note that the comoving observers move outward in the 
wall region in the GR case. However, when a = 10~^°, the SR and GR results are quite 
similar for the entire collapse. In addition, the peak area looks virtually unchanged in 
both cases. This is easy to explain. Using Eq. ( ^4.6| ) and ([4. 81), the peak changes on time 
scales .27 and 5.4 for the GR and SR cases, respectively. And in the peak area, the time 
scales of change for the GR and SR voids are .055 and .18, respectively, using Eqs.(|4.5|) 
and Q. 

In Figure 5 we plot p and F as a function of time for comoving observers in the peak 
area. These voids are identical to those of Figure 3 with a = 10~^. Again, the solid and 
dashed hues show the results in the GR and SR cases, respectively. Note that the energy 
density for the GR void in the peak and peak area has decreased significantly by t ~ 1.1, 
in agreement with the time scale .055 calculated above. In addition, p keeps decreasing 
for all three observers. In the SR case however, the pressure at the peak hardly changes, 
while that in the peak region has decreased substantially by t ~ 1.3, in agreement with 
the time scale .18 calculated above. In addition, in the SR case the energy of the peak, 
c^(T„(ty — l)/2, is substantially less by t ~ 2. In the GR case however, T^(t) starts to 
decrease, then stops when the peak energy density decreases substantially. Therefore, if a 
superhorizon-sized void collapses in more than the outside Hubble time, the wall energy 
probably does not decrease, as it does in the SR case. 

2. Collapse Times for Superhorizon-sized Voids 

We now calculate the collapse times for many different superhorizon-sized voids. 
These times depend on the initial wall energy, c^(r„^ — l)/2, the relative energy density, 
a, and the wall thickness Ai?„. Again, we choose a void which is 25 times the Hubble 
radius, with relative wall thickness AR^/R^ = 1/50. We set AR(t{) = .25, P = 1.7, and 
determine Ur from Eqs.( ^.20D , ( ^.17] ) and ( |3.19D by setting b = 1. We present our results 
in Table 3. The stars following some of the collapse times mark those voids where there 
is significant distortion observed in the peak area just before the collapse. It is interest- 
ing to compare this with the time scales for change in the peak area. Using Eqs. (|4.5| ) 
and ([4.3|) , the time scales for change are given also. As expected, for a given value of 
r^, the peak area is distorted unless the collapse time is smaller than the time scale for 
change in the peak area. In addition, if the peak is distorted and is increased further, 
the collapse time increases. For a = 10~^, a = 10^^ and a = 10^^°, the collapse time 
increases by 21%, 4% and < 1%, respectively, for F„ = 6 to 50. This is because more 
fluid moves out of the peak shock) rather than into the void, weakening the 

shock. (The outgoing shocks can be observed in Figures 6-7). This increases the shock 
formation time, and decreases the relative potential of the void because $in decreases 
slightly when an outgoing shock is present also. This can be seen by rewriting Eq. ( |2.9| ) 
as 

V Pin / \uj + lJo 'yp + Q) 

Q is non-zero only in a shock, and is peaked at the center of it. Divide the shock front 
into two portions separated by the position where Q' = and p = pq Qq. Since p > Q 
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in the shock, for an inbound shock the contribution to the integral is larger when p < pq 
{Q' > 0) than when p > pg {Q' < 0). Therefore when the shock is inbound (outgoing), 
$in increases (decreases) from its isentropic value. As an example, if we assume that 
both inbound and outgoing shocks are strong, using Eqs.( [A.29D and ( |3.10| ) we find that 
the relative potential of the void is 




$i„ = «-^M— — ' (4-12) 



where pb^i {pb,o) are the energy densities behind the inbound and outgoing shocks, respec- 
tively, and a = Pmiti) / PoutiU) ■ Thus, the outgoing shock partially cancels some of the 
void's relative potential. 

In the last two rows, we calculate the theoretical collapse times for each a using Eq. 



( [4.9| ). Note that $in = a is the isentropic result, and $in = .6a is the relativistic 



shock tube result from Eq. ( |5.24| ). Therefore, the numerical collapse times are nearly that 



for the isentropic "no-shock" result. This is because the void wall is too thick, which 
results in the production of weaker shocks, large shock formation times and a smaller 
time dilation effect. 

We can show this by calculating the collapse times as a function of the width Ai?w 
for a superhorizon-sized void with R„ = 50, a = 10~^°, = 6, AR{ti) = A_R„/4 and 
6 = 1. For AR„ = 1, .5 and .25, we set k"^ = 1.7, 2.2 and 3.0, respectively. The results 
are presented in Table 4- In the last column, we calculate the exponent p using Eq. (|4.9| ) 
with the assumption that $in = .6a~^. As AR^ decreases, the collapse time decreases, 
and p increases, as expected. 

As a final remark, because the time scale for change in the peak area of a SR void is 
greater than that for a GR void with the same initial conditions, if a GR void collapses 
in much less than the time scale for change in the peak area, we expect the GR and SR 
voids to look the same and be unchanged in the peak area at the collapse time.[| This is 
observed numerically for the voids listed in Table 3. The fluid which is in motion at the 
first-crossing time has been accelerated into the void by fluid forces. Therefore, it is not 
surprising that the solutions for these GR voids are nearly identical to those of SR voids 
evolved under the same initial conditions. Thus gravitational effects are not important 
at the first-crossing times for superhorizon-sized voids with Atcoiiapsc ^ Atpcak area- 
3. Peak Distortions: Fluid or Gravitational effects 

From the above discussion, it is clear that the wall of a void will be distorted at the 
collapse time if is large enough, even if the collapse occurs in less than an outside 
Hubble time. It is therefore important to find out if this distortion is caused by non- 
gravitational (fluid) or gravitational effects. Thus, we compare the GR solution before 
collapse with that of the SR solution with the same initial conditions. We choose the 

initial conditions such that Atpeak area < AtcoUapse- 

''It is possible to find initial conditions such that the collapse time is between the GR and SR 
time scales for change in the peak area. However, in order for this to happen for collapse times 
much less than the Hubble time, Ai?„/i?w must be extremely small. For example, using Eq. (^^), 
if Atpcak aroa/i?out(^i) - -001 and .005, we require AR„/R„ ~ 2.3 x 10""* and 1.9 x 10"^ xhesc voids 
have extremely thin walls, and will not be considered here. 
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We evolve a GR and a SR void with = 50, AR^ = l,a= 10"^°, r„ = 200 ~ U^/c 
and (Tr = l.SSAi?^. In addition, k"^ = 1.7 and AR{ti) = .25. It is found that the collapse 
time for the SR void is At^ = -18, which is 31% larger than the GR collapse time of .13. 
In Figure 6, we plot p and F as a function of R{t, r) for the GR void (dotted line) and 
SR void (dashed line) at t = 1.06. The initial void configuration is shown as the solid 
line. Note that the GR void's shock is ahead of the SR void's shock (24 < R < 26), 
which explains why the SR collapse time is longer. Since the time scale for change in the 
peak area, .0052, is much smaller than the collapse time, we expect the large distortions 
seen in the peak area (45 < R < 60); the peak has moved outward and is thinner. In 
addition, the peak energy is smaller, which is due partly to the spherical damping effect 
of the outbound wave. However, the important and surprising result is that the GR and 
SR voids look very similar at this time. It is true that the energy in the void wall should 
drive the fluid out regardless of whether gravitational forces are included or not. However 
once the fluid starts expanding outward, one can imagine gravitational effects adding or 
subtracting substantially to this movement. Apparently, this does not happen. 

We would like to know how robust this result is. Therefore we examine voids with 
larger and smaller collapse times and for which peak distortions are large. We also require 
Fw ^ i?w/(v^cif^t(ti)) ~ 18 so that the initial velocity in the wall is the same for both 
GR and SR cases: ^ cF^.Q If we choose smaller values for F^, then it is not possible 
to distinguish the non-gravitational from the gravitational effects in this manner. 

For a void with R^ = 50, AR^ = 1, a = 10"^ F„ = 50, ar = 1.73Ai?„, = 1.7 
and AR{t]) = .25, the GR and SR solutions are very different at t = 1.4. However, a 
substantial fraction of an outside Hubble time has passed. We therefore choose a void 
which collapses much more quickly and for which F^ is very large. 

We choose a void with R^ = 50, Ai?„ = 1, a = lO'^^, F„ = 500, ar = 2.05Ai?^, 
k"^ = 1.7 and AR{t]) = .25. Figure 7 shows that the solutions are nearly identical at 
t = 1.01, even though the peak is fairly distorted. In addition, the collapse times for the 
GR and SR void are Ate = -019 and .021, respectively, which are not very different. Thus, 
the gravitational contribution to the solution is smaller than that for the superhorizon- 
sized void in Figure 6, even though the peak distortions are comparable. This is because 
the collapse time is even smaller. 

In conclusion, if a superhorizon-sized void collapses in much less than an outside 
Hubble time, gravitational forces can be neglected in the equations of motion up until 
the collapse time. This is true even if the fluid in the peak area moves substantially 
during this time. 

V. Lagrangian Shock Tube as an Approximate Solution 
for the Collapse of a Superhorizon-sized Void 

We have found that if a superhorizon-sized void collapses in much less than an outside 
Hubble time, then even if the void wall moves outward significantly, the solution at the 
collapse time looks nearly identical to that for a SR void evolved under similar conditions. 
This is an important result, since it shows that the force moving the fluid out of the peak 

""If Fw ^ Rvj/{V^cH~^^{ti)), then the waU velocity Uv, is much larger in the GR case as compared to 
the SR case, causing the time scale for change to be smaller. 
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area is not due to gravity, but is due to pressure and velocity gradients which are present 
in both cases. (The force moving the fluid into the void is non-gravitational also, for the 
same reason). More importantly, it shows that gravitational forces are not important 
until after the void collapses. This leads us to search for collapse solutions to the special 
relativistic equations (i.e. = 0). 

Consider an uncompensated void with a small wall velocity. As the void starts to 
collapse, the solution is approximately slab (linearly) symmetric. Not until the shock 
nears the origin will spherical effects come into play. If in addition the wall is thin 
(I logio a| A_R„/_Rw ^ 1); then it acts like an initial discontinuity in pressure, energy 
density, etc. The solution for this is called a shock tube. When the solution is assumed 
to be spherically symmetric with analogous initial conditions (i.e. the problem we would 
like to solve), a similarity solution does not exist. This is because another length scale, 
the void radius, R^, has been introduced into the problemj^. Thus, the slab shock tube 



solution derived here will only be important for uncompensated voids with thin walls, 
small wall velocities, and for which AtcoUapsc ^ H~J^(ti). 

The general relativistic equations used here do admit similarity solutions ||12|. We will 
not solve these equations here however, because we cannot find similarity solutions for our 
problem. In Appendix A, we review the conditions satisfied over a contact discontinuity 
and shock for a general relativistic fluid0. In section VA, we derive the slab similar- 
ity solution for relativistic fluids in this Lagrangian coordinate system. We then attach 
this solution to the contact discontinuity and shock solutions, which are also similarity 
solutions, when the fluid in front of the shock is relativistic and nonrelativistic. This 
corresponds to the fluid in the "void" being relativistic and nonrelativistic, respectively. 
The relativistic (shock tube) solution was derived in a Eulerian coordinate system, which 
can be obtained from the Lagrangian description by a spacetime coordinate transforma- 
tion|T^. Finally, we compare these solutions to the numerical solutions obtained from 



the collapse of relativistic or nonrelativistic voids. 

A. Relativistic Similarity Solution 

In this section, we derive the similarity solution for a relativistic fluid (with no viscos- 
ity) in the slab symmetric limit. This limit corresponds to the evolution of a spherically 
symmetric fluid distribution which is "far" from the origin. Mathematically, we require 
\U/R\ < \U'/R'\ to be satisfied. 

We define X{t,x) = R{t,r) to be the location of a comoving shell with label x = r. 
In order that the solution be a similarity solution, all functions can only depend on the 
variable 

z^^, (5.1) 

where t = is the initial time. Thus, p = p{z), p = p{z), U = U{z), $ = $(2;) and 
xi^) = X/x. {X = R and x = r in the slab limit). We can then rewrite Eqs.( |2.4|) - (|2.5|) , 
(g§ and ( CT ) as 

^dp ^ Hp^ + p) dU 
dz dx/dz dz 
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dz ' ' 

^dU_ _ c^r^$ dp 

dz {p + p)dx/dz dz 

^ ^ (5.5) 

dz p + p dz 

dn $n dU 

dz dx/dz dz ' 

where x = = ^ — p/n — c^. We now assume that e/c^ ^ 1 and that the 



equation of state is p = ujp. Eqs.( p^) - (|5.6| ) can then be solved to obtain 

$f/ = ^cs^T + x (5.7) 
p/po = (^Ao)^'^"^/^'""^ = (toA)^'+^^/^'-"^ (5.8) 

Ulc - - ^ (5 9) 

/ " 2(p/po)/^ ^^-^^ 

$ = $o(po/p)"/('+") (5.10) 
n = n,{p/p,f'^^+-^ (5.11) 

and e = p/n, where /5 = ±yZJ/(l + cu), C5 = y/uc is the speed of sound, and where the 
subscript "0" denotes the value a shell has upon entering the region where this solution 
holds. Eqs.( p.7D -( p.ll| ) describe a rarefaction wave going to the right (left) when the 
upper (lower) sign is chosen. In addition, the position of a shell is given by 



t 



±2cs-{c±cs){l-{p/po?^} 



2(p/po) 



/3 



(5.12) 



Therefore, we have solved for all functions in the rarefaction wave in terms of (p/po). 

B. Shock Tube Solution 

In this section, we find the solution for a special relativistic fluid with an initial 
discontinuity in the pressure. The side with the larger pressure is relativistic, and the 
fluid on the other side is relativistic or nonrelativistic. These solutions will be used to 
model the collapse of super horizon- sized voids in the early universe. 

Suppose at t = ti, the fluid is at rest everywhere, with an initial discontinuity at 
position X„(ti). On either side of the discontinuity, p, n, p and e are constant. The 
subscripts 1 and 2 represent the fluid quantities for X(ti) < X^(ti) and X(ti) > X^(ti), 
respectively. We choose the coordinate system where P2 > Pi, ^2 > ei and n2 > rii. The 
equation of state is p = up in region 2, and p = up {p = 0) in region 1 if the fluid there is 
relativistic (non-relativistic). For t > ti, a shock will be formed which travels into region 
1 (undisturbed fluid), a rarefaction wave will be formed which moves into region 2 (also 
undisturbed), and a contact discontinuity will form in between. There will be 5 distinct 
regions in all. The shock is located at Xd(t, Xgh), region 3 contains the fluid behind the 
shock, the contact discontinuity is located at X^it, x), and the rarefaction wave is located 
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in region 4 for Xb{t,x) < X < Xa{t,x). Figure 5 is a sketch of pressure versus distance 
for a relativistic shock tube, with the regions and boundaries labeled. 

We will assume that the fluid in the rarefaction wave is relativistic regardless of the 
equation of state in region 1, so that es'/c^ ^ 1. Doing this allows us to use the similarity 
solution obtained in VA for the rarefaction wave. And because > e^, as we will see in 
a moment, it therefore follows that es/c^ ^ 1 in both cases. Then, the solution in region 
3 is given by Eqs. (|A.2CI|) - (|A.24 ) if region 1 is nonrelativistic, and Eqs.( A.26 )-( A.30 ) if 
it is relativistic. In either case, the upper signs are used because the shock moves in 
the direction of decreasing X. Note that these are all given in terms of ps = pb- The 
conditions across the contact discontinuity give $3 = $3/, ^3 = py, = py and U3 = Uy 
from Eq. ( |A.15| ). Finally, the rarefaction wave solution is given by Eqs.( |5.8|) - (|5.11|) with 
the upper signs and p2 = Po- 

1. Relativistic, Thin- Walled Uncompensated Void 

We first examine the case where the fluid in region 1 is relativistic. Since Uy/c 



\2(5 



l)/[2{p^/ P2Y] from Eq. (|]9D, we can use this and Eq. (|X]2|) to find the 



equation relating p = ps/ pi, (3 = ^/to/{l + to) and 5 = (pi/p2)^: 
25/5(p-l)/-i/2 = l-5y^. 



(5.13) 



An iterative method is then used to find p. The solution in region 3, 3' and 4 is given 
by Eqs.([A:2^)-(|A:30|), Eqs.(p:9D-(pTTD with p = P3 and po = P2, and Eqs.(|5:8|)-([5Tl|), 
respectively. 

We first calculate the location of the boundaries between the five regions: Xa, Xb, Xc 
and Xd- Set $2 = 1- We can calculate Xa by noting that the velocity is zero there and 
by using Eq. (|5.7|). Xb is found by using Eq. ( ^.7| ) and ( [A. 26] ). Since U{z) and $(z) are 
constant across the contact discontinuity and in regions 3 and 3', we can integrate Eq. 
( ^.3|) to obtain the position of a shell in region 3 or 3': 



X{t,x) = ^Uit-to)+Xo. 



(5.14) 



Here Xq is the location of a shell when it enters region 3 or 3' at time to. Eq. {^J§ 
also applies to the location of the shock and contact discontinuity, with Xq = X^(ti) and 
= t\. The locations of the boundaries are then 



Xa{t) 

Xb{t) 

Xc{t) 
Xd{t) 



X^{t{)+cs{t-U) 
X„(ti) + (t-ti)$3/5- 



^p 



-C + Cs\ 1 + 



(5^{-p-iy 



Xw(ti) + (t-ti)$3f/3 

X„(ti) + {t- t,)^^U,^ = X^it,) + (t - ti)(l + ^] 



p-1 



-$3t/3. 



(5.15) 

(5.16) 

(5.17) 
(5.18) 



where the shock speed is ^zUsh = -Rsh^sh = -^sh^sh from Eq. ( |A.30|) . Note that the 
rarefaction wave moves out at the speed of sound. In addition, from the discussion 
following Eq. ( |A.30| ), if the shock is strong (i.e. p 3> 1), the shock moves into region 1 at 
the speed of light: Xd{t) = X^{ti) - c(t - ti)$i. 
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Now we trace the location of all Lagrange points in the shock tube in terms of z = 
x/{t — A shell with coordinate x and position X initially in region 2 [X > X^{ti)) 
will enter the rarefaction wave at time ta = ti + {X — X^{t\)) / cs- Using Eq. ( p.l2[ ) and 
setting Za = x/ (ta — ti), its position as a function of time in that wave for z < Za will be 



Xi{t,x) =X„{ti) + {t-ti) - 



cs 



cs + c 




(5.19) 



Here, the subscript "4" shows that the shell is in region 4. The energy density of this 
shell decreases until it reaches pa, at which point it enters region 3'. This occurs at time 
h = ti + {ta — ti){p2/ p-iY^''^^^^^'^'^^ ■ Its position at this point is Xh(tb,x) and is given by 
Eq. or Eq. ( [5. 191 ) with t = ti, and z = Zb = x / {tb — ti) . Using Eq. (|5.7D, its position 

in region 3' for z < Zb is 



Xsit, X) = X^iti) + (t - ti) I $3t/3 - - ( $3f/3 



— ^witi) 
tb — ti 



(5.20) 



This Lagrange shell can never cross into region 3 because the contact discontinuity moves 
with the fluid. 

Finally, if a shell initially has the position X < X^{ti), then it will enter region 3 
when the shock reaches it. This occurs at time td = ti + {X — X^{ti))/ {^2,Ush)- For t > td, 
the location of this shell in region 3 for \z\ < \zd\ is 



X,{t, X) = X^ti) + {t- ti) $3f/3 - - $3f/3 



Zd 



td — ti 



(5.21) 



where Zd = x/{td - ti). 

In the strong shock limit (p ^ 1), p is given by an analytic expression. The solution 
then is 



P3 

Us 
c 

$1 

$3 

£3' 
£3 



1+UJ 

— w 



1.2pi«-^64 
-.47(5-232 



2^ 



' (l+v^)2 



.96$2a" 



.134 



e2($2/$3) 1.01e2a-i3^ 
ei^/ujJp3/pi .62eia"-23^ 



(5.22) 
(5.23) 

(5.24) 

(5.25) 
(5.26) 

(5.27) 



where a = pi/p2 and a = [(1 + u) / {2y/uj)Y^^^^^^'^ , and where the arrow evaluates 
the expression when u = 1/3. We note that in general €3 7^ £3/. This is because 
energy is dissipated over the shock, increasing the specific energy for shells behind the 
shock. Suppose we take our initial state to be an isentrope so that ei/e2 = {pi/p2)''^^- 
Then, £3/63/ = .6(p2/pi)'^^^ > 1. Thus the specific energy and number density will be 
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discontinuous across the contact discontinuity, and will diverge as pi ^ 0. In order 
that this solution be consistent, pi/p2 cannot be too small from Eq. ( |5.26| ) because we 
assumed that e^i/c^ S> 1. In addition, since the temperature in region in region 3 is 
T3 = w/xes which we require to be less than the Planck temperature, the energy density 
inside the void cannot be too small. However, we see no reason to restrict a for more 
general solutions. 

The above solution approximates that of an uncompensated super horizon- sized void 
in a radiation-dominated universe. We can calculate the approximate energy density in 
the "void" at times much less than a Hubble time before collapse. This solution gives 
the theoretical maximum limit for the shock tube strength. Of course, spherical effects 
will cause the shock strength to increase. However, our numerical simulations suggest 
that this increase is proportional to the initial shock strength. The ratio of the energy 
density behind the shock to that in region 2 is 

\ 1/(2/3+1) 



P3 
P2 



1 



a 



1.2a 



;r,.536 



(5.28) 



where the arrow evaluates the expression for = 1/3. In addition, the collapse time is 
AtcoUapse — — 2_R„a-^''. Thus as the "void" empties {pi/p2 — ^ 0), the collapse time 

decreases rapidly and the fraction filled before the "void" collides at the origin becomes 
very small, unless spherical effects are extraordinarily important. 

2. Nonrelativistic, Thin- Walled Uncompensated Void 
Now we calculate the solution for a shock tube when region 1 is nonrelativistic: 
Pi = riic^ and ei/c^ <C 1. We assume that the fluid behind the shock and in the 
rarefaction wave are relativistic, however. We use Eq s . (|A . 2 0|) - ([A . 24|) with 63 /c^ ^ 1, 
as well as the same rarefaction and contact discontinuity solutions as was used for the 
relativistic case in the previous section. In particular, we use Eq. ( p.9|) with p = ps and 
Po = P2 in order to obtain f/3/. Since ps = n^es = {{u + l)/uj) piies/c^)^ and f/3 = — ea/c, 
U3 = U3' becomes 



,^3 



UJ + 1 Pi 



UJ P2 \C 



1 



UJ 



Ipl 



UJ P2 VC 



2/3 



(5.29) 



This equation can be solved iteratively to find e^/c^, from which point the entire solution 
can be obtained. We assume now that P3/P2 -C 1. In this limit, we obtain the following 
solution: 

-1/(1+2/3) 



P3 

$1 

$3 
£3' 



{2[a{uj+l)/ujfy 
[{uj + l)/uj] pi el/c' 

<l>i/[(^+l)e3/c2_ - 
€2(^2/^3) 1.0e2a 



.50a 



;t--232 



pia-^'' 



a^+" 



v^(l + V^)2 



(2/3)$2a 



--.366 



<l>2a 

-.134 



--.134 



(5.30) 

(5.31) 
(5.32) 

(5.33) 

(5.34) 
(5.35) 
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where a = pxj P2- This solution depends on a in the same way as it does in the relativistic 
case. The prefactors however, are somewhat different. Using Eq. ( |A.25| ), ( |5.31| ), ( p.33| ) 



and ( |4.9| ), we find a rough lower bound for the collapse time 



where /i is the particle mass. (In general, because Ty < Mp^, we can estimate the lower 
bound on the collapse time using Eq. ( |3.12|) to be AtcoUapsc > c^^Rw p/MpLj- We can 



also trace out the boundaries between regions, and the radii of all comoving shells. Xa 
and Xc are still given by Eqs. (|5.15|) and (p.l7|) , respectively. However, now 

X, = X„(tO + (t-ti)$3f/3(l-v^) (5.37) 
X, = X„(tO + $3f/sh(t-ti)=X^(ti)-c(t-ti)$i. (5.38) 

The rarefaction wave moves out at the speed of sound, and the shock moves in at roughly 
the speed of light. 

C. Comparison of the Collapse with the 
Shock Tube Solutions — Numerical Results 

In this section, we find the numerical solutions to the special relativistic equations of 
motion with spherical symmetry and compare them to the slab-symmetric shock tube 
solutions derived in the last section. For Figures 9-12 and Table 5, we set C = .3, c = 1, 
= 1/3, ti = 1, Gn = 0, U{ti, R) =0 and 47rpout(ti) = 3/8. 

In Figure 9 we show the collapse of a relativistic, uncompensated void with R^ = 20, 
AR^ = .005, AR{ti) = .0025, a = lO'^, eoutiti)/c^ = 10^ and = 2.5. The triangles 
show the pressure and specific energy as a function of R at ti and t = 1.5. As usual, 
the void collapses via a shock moving toward the origin at the speed of light. The shock 
tube solution is plotted as dashed lines. The agreement is excellent, because the shock 
has not traveled very far and \ log^g C(\AR^/R^ <C 1. 

Due to time dilation in the void, light travels a farther distance inside the void than 
outside it in the same amount of coordinate time t — ti. Therefore, the amount of time 
taken to form the shock is very important for relativistic fluids. If the wall is too thick, the 
shock will take longer to form. In addition, the shock strength will be smaller, decreasing 
the void's relative potential and therefore the time dilation effect. This results in the 
solution not equaling the shock tube solution at any time. In Table 5, we show the 
location of the shock, Rsh, at t = 1.5 for voids with the same initial conditions as that in 
Figure 9 but for varying wall thicknesses and AR{ti) = AR^/2. In addition, the value of 
$in is shown. Note that $in(ti) = 3.16 for all voids. It is clear that for this value of a, the 
wall must be thinner than 2AR^/R^ ^ 2 x 10~^ in order that the solution approximate 
the shock tube well. 

In Figure 10 we plot p, e and f/ as a function of time for the void in Figure 9 but 
at t = 5. Note the large distortion in pressure, specific energy and velocity behind the 
shock. This is due to spherical geometrical effects; fluid behind the shock is forced to 
occupy a smaller volume, thereby increasing the number density. This strengthens the 
shock, which in turn dissipates more energy so that the specific energy, energy density 
and velocity increase. In addition, it is seen that the shock is far ahead of the shock-tube 
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shock. This is a consequence of the increased energy density behind the shock, p}j. From 
Eqs.( |A.29| ) and ( |3.10| ), the value of $ inside the void is 

Therefore, as the energy density behind the shock increases, the shock moves further per 
unit coordinate time because (i-Rsh = c^dt. 

In Figure 11(a), we show the shock tube solution for a = 10^^ and X^(ti) = 10. The 
dotted lines are the position of Xa, Xb, Xc and Xd, and the solid lines are the locations of 
comoving observers with X{ti,x) = 11.0, 10.1, 10.0, 8.0, 4.0 and 2.0. Note that the shock 
moves very far per coordinate time interval, due to the time dilation effect. Note also 
that after entering regions 3 or 3', all shells move with the same velocity. In Figure 11(b), 
we show the comoving shell positions again as solid lines. In addition, the dashed lines 
are the positions of comoving shells with the same initial locations, R{t[,r) = X{ti,x), 
in an evolving void. For this void, = 10, AR^ = .1, AR{ti) = .025, a = 10~^, 
^ont{ti)/c^ = 10^ and k'^ = 1.7. Although the void wall is not very thin, the solution 
approximates the shock tube solution fairly well in region 3. However, it does not do 
as well in regions 3' and 4, due to effects stemming from the thick wall. Note that 
the numerical shock collides at the origin at t = 2.2, which is why the velocity of the 
comoving shell with R{ti, r) = 2.01 changes after this time. 

In Figure 12, we plot p, e n and U for the collapse of a nonrelativistic void in a 
relativistic background: ein(ti)/c^ = .38, and eout(^i)/c^ = 5.0. In addition, R^ = 20, 
AR^ = .0025, AR{ti) = AR^/A, a = lO'"^ and A;^ = 3.0. The triangles are the numerical 
solution at t = 1.06.| The general features are similar to the relativistic void solution: 
the void collapses at the speed of light as the wall fluid shocks inward. In addition, a 
contact discontinuity (where p and U are continuous, while e and n are discontinuous) 
follows behind the shock. Finally, a rarefaction wave moves into the region outside the 
void. Note that the shock and contact discontinuity are located at R = 19.33 and 19.47, 
respectively, while the rarefaction wave is located between R = 19.7 and 20.04. 

The dashed line is a plot of the nonrelativistic/relativistic shock tube solution derived 
in section B2. It is clear that this solution agrees well with the numerical solution in 
the outer part of the rarefaction wave (large R). In addition, p, e n and U are predicted 
fairly well behind the shock. (For example, ps/pi = 61 and 63 = 4.2 from the numerical 
solution, while ps/pi = 72 and e/c^ = 4.2 from Eqs. (|5.31|) and (|5.30| )). However, because 



the fluid in regions 3' and part of 4 are barely relativistic (e/c^ ~ 1), the assumptions 
used to derive the shock tube solution are not satisfied in this case. Therefore, it is no 
surprise that the shock tube solution does not match the numerical solution very well. 
At t = 1.06, $in = 11.2 from the numerical solution, while $1 = 19.4 from Eq. ( |5.33| ). 
Note that $in increases from its initial value $in(^i) — Cout/c^ = 5.1 because energy is 
dissipated over the shock. This results in a smaller collapse time. 



"Note that the numerical simulation covers the range R{ti,r) g [19.2,20.1]. In addition, Ai?w is 
chosen to be small enough that the wall acts like an initial discontinuity. 
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Thus, the relativistic shock tube solution approximates the numerical solution of a 

thin-walled void very well for a significant amount of time during collapse. The nonrela- 
tivistic/relativistic shock tube solution does not do as well because the fluid between the 
contact discontinuity and the rarefaction wave tends to be barely relativistic or slightly 
nonrelativistic. However, it is qualitatively (and to some degree quantitatively) correct. 

VI. Discussion 

In this paper, we studied in detail the non-linear collapse of a superhorizon-sized void 
embedded in a background radiation-dominated FRW universe. We find in general that 
a relativistic or nonrelativistic void collapses via an inbound shock at the speed of light. 
This occurs because the pressure and velocity gradients in the lower part of the void 
wall are enormous and therefore accelerate fluid into the void. Because of a time dilation 
effect, this occurs very quickly, and in particular can occur in less than an outside Hubble 
time if the relative energy density inside a void is small enough. This is true even when 
the outward velocity of the wall is enormous. A large wall velocity does drive some of the 
fluid in the wall out behind a shock, so that less fluid ends up in the void. This causes 
the strength of the inbound shock to decrease, and lessens the relative potential of the 
void. Thus the collapse time increases somewhat. However, the void can still collapse in 
less than an outside Hubble time. This is important, since it emphasizes the point that 
seemingly nothing stops the void from quickly collapsing in on itself. 

As the void wall gets thinner, the shock strength increases and therefore the collapse 
time decreases. This behavior does not extrapolate to inflnitely thin walls, however. 
When the wall becomes sufficiently thin so that it acts like an initial discontinuity, a 
maximum value for the shock strength (and therefore a minimum value for the collapse 
time) is reached. Up to the point when the void has partially collapsed, the solution 
can be approximated by the shock tube solution. This holds for an uncompensated, 
superhorizon-sized void which collapses in much less than an outside Hubble time and 
for which the fluid in the void is relativistic or nonrelativistic. If the fluid in the void 
is nonrelativistic, then then the collapse takes longer to occur for the same value of 

Pin{ti)/Pout{ti). 

If the wall velocity is small and the void collapses in less than an outside Hubble time, 
then the solution at the collapse time is virtually the same as that when gravitational 
effects are neglected. This is because fluid in the lower part of the wall is accelerated 
into the void by pressure and velocity gradient forces, not gravitational forces, and the 
fluid conflguration in the peak area is virtually unchanged. 

When the wall velocity is large however, it is not obvious a priori whether or not 
gravitational forces are important in moving the peak area fluid out of the void wall. In 
this case, fluid in the peak area of the wall moves out in less than an outside Hubble 
time. In addition, the time scale for this to occur is proportional to AR^^^^ /Ty,, where 
A.R^ is the wall thickness, and is approximately the "gamma-factor" or "velocity" 
of the wall when > / H~^^{ti) . (II — 0, F^ is the relativistic gamma-factor 

of the fluid, and equals the "velocity" when r„ ^ 1). Therefore, if r„ is large enough, 
fluid in the peak area will shock outward during the collapse, distoring the shape of the 
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wall in the process. It is found however, that in general this does not change the collapse 
time by very much, given a fixed value of p-miti) / Pontit\) ■ In addition, it is found that 
gravity is unimportant in shaping the distortion of the fiuid configuration in the peak. 
Therefore, gravitational contributions to the solution at the collapse time are negligible. 
This means that only non-gravitational forces are driving the fluid into and out of the 
void wall, even though the void is super horizon- sized. (Although it is clear that pressure 
and velocity gradient forces must contribute to driving the fluid out of the peak area, one 
might envision gravitational forces changing the conflguration by substantially redshifting 
the wall pressure, as it does when the collapse occurs in greater than an outside Hubble 
time. Apparently, gravitational forces do not act quickly enough to do this.) This sounds 
surprising, especially since superhorizon-sized voids are larger than the Hubble radius and 
therefore must be influenced heavily by gravitational forces. It is a question of time scales, 
however. Because time is dilated in the void, if the collapse occurs in much less than an 
outside Hubble time, gravity has not had a chance to substantially change the solution 
(pressure distribution, shell locations, etc) at the collapse time in this case. Gravitational 
forces should be extremely important after the void collapses, of course, which was not 
studied in this paper. Then, they will try to pull fluid lumps together, while fluid forces 
will try to pull fluid lumps apart. 

There are many questions still left to answer. What happens to a superhorizon- 
sized void after it collapses? Now that we have a code which implements numerical 
regularization to prevent instabilities near the origin, it will be possible to explore this 
late stage numerically. In SV, we found that after a special relativistic void collapses, 
the pressure near the origin is much greater than the pressure outside the void, creating 
an overdensity. Because this lump drops off with distance from the origin, this large 
pressure gradient accelerates fluid away from the origin. What happens if gravity is 
added to this problem? It appears likely that a black hole might form if the wall velocity 
is not very large. If it is very large, enough fluid in the peak area might be pushed out 
of the void initially so that this is prevented. In this case, a mass-energy deflcit will 
occur in the void initially. If this mass-energy continues to be pushed away from the 
void, thermalization and homogenization of this void might take a long time to occur. 
However, as this expelled fluid slows down, its mass-energy M will decrease. In addition, 
at the collapse time the velocity in the former void region is negative, so its mass-energy 
will increase. Thus, it will be important to analyze the long-term solution to see if there 
is a mass deflcit created where the void was located. In addition, density waves will be 
created, and must also be followed for their long-term behavior — do they become growing 
density perturbations, or remain outgoing local waves? A companion paper will explore 
these questions 0]. 
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Figure Captions 

Fig. 1: Relative error in energy density, (p — Phom)/Phom, where phom is the FRW 

solution, versus RJti/t for 3 different values of /. The improved code does much better 
than that used in SV. 

Fig. 2: The relative difference in pressure as a function of initial radii for the collapse 
of GR and SR voids with nearly identical initial conditions. = 1 for both. The closed 
triangles, open squares, and dashed line show the results for a — 10~^, a — 10"^ and 
a = 10"^*^. As the collapse time becomes much less than an outside Hubble radius 
{a < 10~^°), the GR and SR solutions are very similar. 

Fig. 3: Same as Figure 2, but for F^ = 6. 

Fig. 4: The evolution of six comoving observers as a function of time for the voids 
from Figure 5with a = lO-^(top) and a = 10-^°(bottom): R{ti) = 30, 35, 40, 45, 48, 50, 
55 and 60. The solid (dashed) hues are the GR (SR) solutions. These solutions arc very 
similar for a — 10~^° because the collapse occurs in much less than an outside Hubble 
time. 

Fig. 5: 47rp(t) and F(t) for 3 comoving shells in the peak area, R{ti) = 50, 49 and 48 
(top to bottom) for the voids from Figure 3 with a = 10~^. The GR (SR) void results 
are shown as solid (dashed) lines. Note that the peak energy density redshifts away quite 
quickly for the GR, but not for the SR, void. 

Fig. 6: 47rp and F as a function of radius for a superhorizon-sized void 25 times the 
Hubble radius and with a — 10~^°. The initial distributions are shown as solid lines. The 
numerical GR (SR) results are plotted as dotted (dashed) hnes. Note that the peak has 
moved outward substantially at At / H~^^{ti) — .03 for both voids, but that the solutions 
are very similar anyway. 

Fig. 7: Same as Figure 6hut for voids with relative energy density inside the void 1000 
times smaller. Note that the GR and SR solutions are nearly identical at At/H~^^(ti) — 
.005, even though the peak area has shocked outward substantially. 

Fig. 8: Schematic diagram of the shock tube showing the different regions. The 
shock is located at Xd, and is moving into region 1, and the outer edge of the rarefaction 
wave is located at Xa, and is moving into region 2. 
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Fig. 9: Pressure and specific energy as a function of radius for a thin- walled SR 
void. The initial distributions are the solid lines, the numerical solution is the triangles, 
and the shock tube solution is the dashed lines. Note that the numerical solution is well 
approximated by the shock tube solution. 

Fig. 10: Same void as in Figure 9, but at a much later time. The numerical solution 
is given by the solid lines in this figure. The pressure, specific energy and velocity behind 
the shock is much larger than the shock tube solution. This is due to spherical geometrical 
effects. The upper part of the rarefaction wave however, is unaffected by these effects. 

Fig. 11: (a) Location versus time for the shock tube solution with a = 10~^ and 
— 10. The position of the 4 boundary locations are the dotted lines, and the solid 
lines are the locations of comoving points with X(ti) = 2, 4, 8, 10, 10.1 and 11. (b) 
Shock tube solution (solid lines) and numerical solution (dashed lines) for the same 
comoving points from Figure 11(a). The two solutions coincide well in region 3, but 
not as well in the wall area {X{ti) = R{t{) = 10, 10.1) because the thick wall evolves 
differently than that given by the simple similarity solution. 

Fig. 12: Solution for an uncompensated, nonrelativistic void in a relativistic back- 
ground (triangles). The dashed hues show the nonrelativistic/relativistic shock tube 
solution. 

Fig. 13: Percent change in e, F and n as a function of shock strength for a 
general relativistic fluid. Note that e is underestimated, while the other quantities are 
overestimated, at the shock. 



30 



Appendix A: Contact Discontinuity and Shock Jump Conditions 

In this section, we briefly re-derive the jump conditions for general relativistic shocks. 
This closely follows Ref . . In addition, we extend the derivation to include the conditions 
satisfled over a contact discontinuity. Let the variables a and h represent labels for 
comoving observers ahead and behind a discontinuity, respectively.^ If each observer 
measures the invariant interval separating the same two events on the world surface of 
this discontinuity, using Eq. (|2.1|) , we obtain 



[ds'] = l-c'^'dt' + A'dr' + R'dn'] = 0, (A.l) 

where [G] = Ga ~ G^. Because the discontinuity is radial, we can choose the two events 
to have dr = dt = but dQ 7^ 0. Therefore, R is continuous across the discontinuity: 
[R] = 0. Now consider two events with dr and dt ^ 0. Then from Eq. ( |A.1| ), 

[c^*^ - Ml/n^] = 0, (A.2) 

where Vgh is the position of the discontinuity in comoving coordinates, drgh/dt is the 
discontinuity "speed", and Msh = f {dr^J dt)/ {AtxR'^) from Eqs.(|2]|) and (IJ). If the 
discontinuity is a contact discontinuity, then rgh = const, or dr^^/dt = 0. 

It can be shown that the conditions satisfled over the discontinuity for the Schwarzschild 
metric (but not for the comoving metric) are [n], 13] 



m^dg/dx"] = 0, (A.3) 

where g is the equation for the world surface of the discontinuity. Therefore we flrst flnd 
these expressions in the Schwarzschild frame and then transform back to the comoving 
frame. The Schwarzschild metric is 

ds^ = -c^A^dT^ + B^dR^ + R^dn^, (A.4) 

where R is now the "Eulerian" coordinate radius. Since i? is a coordinate, [R] = over 
the discontinuity, and using the same argument as above for the two observers, 

[A^ - B^S^] = 0, (A.5) 

where S is the discontinuity "speed" in this frame: cS = dRg^/dT. For the Schwarzschild 
metric, the solution for the metrics functions are well known: B~^ = 1—2MG/ {Rc^) = A? 
where M{R,T) = Attq-^ pR^dR. As long as p is not inflnite in the discontinuity, the 
mass is continuous across it, [M] = 0. Then [B] = [A] = 0, and from Eq. ( |A.5| ), [S] = 0. 

We now derive these conditions in the Schwarzschild frame. If a discontinuity is 
located at position _Rsh(^) at time T, the equation for its world surface is g = RshiT) — 
R = 0. In addition, the perfect fluid stress-energy tensor in the Schwarzschild frame^ is 

= c-\p + p)g',^u'^u'^ + pg'^J'^^ = nWg'^,u"^u" + pg'^,g"^\ (A.6) 



*Here, "discontinuity" refers to either a contact discontinuity or a shock. 
"We denote the metric and 4— velocity in the Schwarzschild frame by primes. 
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and the conservation of mass equation is [nu"^dg/dx^] 
( |A.3| ), the junction conditions become 



0|T1],[P1. Using this and Eq. 



-A^u^'fnW + p)S + nWA^'^u"" 



nSu 



IT 



nu 



iR 





0. 



(A.7) 
(A.8) 
(A.9) 



We would hke to express these equations in terms of the comoving metric functions. 
First, we can rewrite the speed in the Schwarzschild frame as cS = {Rt + Rrrsh)/(Tt + 
TrTsh)- In addition, we can relate the metric functions from the comoving frame to 
those in the Schwarzschild frame by g''^" = {dx'^ / dx") {dx'" / dx^) g"^. We then obtain 

-'^TtRt = A-^TrRr. In 



{cA) 



-2 



(C$) 



A-^T' B 



-2 



-(C$) 



-^Rj + K' 



"^Rl and (c*)- 



addition, the "mass-energy" M(r, t) is continuous across the discontinuity, since [M] = 
ATic-'^[pR^]dR 0. With this and Eq. pAl| ), we find that [T^ - U^/c^] = 0. 

We can use these relations to calculate the 4- velocity as measured in the Schwarzschild 
frame. Since the 4-velocity in the comoving frame is = (— c$~^, 0, 0, 0), and the ve- 
locity transforms as u''*' = {dx''^ / dx'^) u" , u'^ = — $^^(cT, R, 0, 0), where T = dT/dt and 
R = dR/dt. Using the normalization condition w 



rX 



— C 



c'AHu'^Y + BHu 



iR\2 



along with U = -R/$, we can rewrite the 4-velocity as u''^ = — (A^^ v5^f7^ + c?, U, 0, 0). 
In the special relativistic limit {Gn = 0), A = B = 1. The 4— velocity of a fluid particle 
is then 



u' 



rX 



-((t/2 + c2)V2,f/,o,o) = -(cr,rt;,o,o), 



(A.IO) 



where we have defined v = U/T so that F = 1/ y^l — (f/c)^. In the special relativistic 
limit, then, v is the fluid particle's radial velocity, F is usual gamma-factor (i.e. energy 
per particle mass), and U is the particle momentum per particle mass. 

Across a contact discontinuity, which is a discontinuity with = 0, Eqs.( [A.2D , ( |A.7| )- 

and [S] = 0, 



[$] 

\pU/T] 

[P] 
[U/T] 










(A.ll) 
(A.12) 
(A.13) 
(A. 14) 



respectively.^ Thus the pressure is continuous across a contact discontinuity, as it is 
in the nonrelativistic case. Using Eq. ( |2.11| ) and ( A. 14 ), we find that {U^ — U^){1 — 
2GnM/{Rc^)) = 0. If 2G'ArM/(i?c2) ^ 1, [U] = 0, or the tangential velocity is continuous 
across the contact discontinuity, as it is in the nonrelativistic case. If 2GnM/{Rc^) = 1, 
however, the entire formalism breaks down because then B = oo. This occurs at the 
Schwarzschild radius. Thus, across a contact discontinuity. 



[$] = 0, [p] = 0, [U] = 0. 



(A.15) 



^Eq. (jA^) gives [0] = 0. 
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We now find the jump conditions at a shock where rsh 7^ 0. We rewrite Eqs.( |A.7] )- 
( |A.9| ) and [5*] = in terms of the comoving metric functions using the fact that T^/Tt = 
c-2f/A/(r$) and [fsj = 0. We obtain 



[t/Msh/(nc2) + $r] 
[c2Mshr(l + e/c') - 
[Mshf/(l + e/c=^)-pr$] 






0. 



(A.16) 
(A.17) 
(A.18) 
(A.19) 



Eq. ( |A.2| ) and Eqs . ( |A . 1"6|) - ( |A . 1 9|) make up the required shock conditions. One of them 
however, is redundant. It can be shown that in the nonrelativistic hmit, the above 



conditions reduce to the Lagrangian shock jump conditions given in Ref.|15 



We are interested in the case where the fluid in front of the shock is either nonrel- 
ativistic or relativistic. First, consider the nonrelativistic case, which was considered 
previously^. We set Gn = and ea = Pa = Ua = and take pb = ojEbUb- Using Eq. 
(^) and Eqs.( [A:TBO -( |AlUD , we find that 



Ub/c 



Ub 



b 



2 + uJ+{uJ + l)eb/c^] /uj 

^ n 



l + {uj + l)eb/c2 
1 + tb/c^ 



sh' sh 



eb/c^ f2 + uj+{u + l)eb/c^ 



\ eb/c^ + 2\ 1 + (1 + uj)eb/c^ 



(A.20) 
(A.21) 

(A.22) 
(A.23) 

(A.24) 



Thus, all quantities behind the shock can be expressed in terms of the specific internal 
energy behind the shock, eb- Note that when the fluid behind the shock is relativistic, 
Cfe/c^ ^ 1, the shock moves at the speed of light, -Rgh^sh = =Fc$a- 

The energy density ahead of the shock cannot be zero since > 0. If the mass 
density were zero, then = 00. But this is impossible since ueb = Tb/^, where /x is the 
particle mass, leading to an infinite temperature behind the shock. Therefore, we require 
Tb < MpL in order that the classical field equations be valid. Then, if e^/c^ ^ 1, the 
energy density in front of the shock is bounded from below: 



Pa = riaC^ > [CUV(CJ + 1)] Pbifl/MpLf 



(A.25) 
leads to 



In addition, it can be shown using Eqs.( |A.17|) and ( |A.18|) that setting 
f/^ — = 0, in direct contradiction with Eq. ( A.20 ). 

When the fluid in front of the shock is relativistic, we end up with a completely 
different set of equations. We first set Gat = and Ua = 0. In addition, p = up and 
p = ne. Then Eq. (^) and Eqns (|A:TBD -( [A:T^ yield 



Ub 



b 



Pb - Pa 
+ 1 y/PaPb 



(A.26) 
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V 



{ujpa + Pb){pa + Ujpb) 
\ (W + lYpaPb 

Ub _ (1 + Uj)TbPb 
Ua Ujpb + Pa 



(cjpa + Pb)pa 
Pb 



Pb - Pa 



-Wb. 



(A.27) 
(A.28) 

(A.29) 
(A.30) 



0. 



We note from Eqn ( [A. 291 ) that 

In the strong shock hmit, when pb 3> Pa-, the shock moves at the speed of hght: 
R'sYJ'sh - ±c$a- In addition, Ub/ua ^ ^ pb/{ujpa) and efe/e^ ~ ^upb/pa- Thus, n and e 

do not scale as p'^/(i+'^) as they do when Q = They always scale as Pbl Pa-, regardless 
of the value of lo. 

In order to see how accurate our code is at the shock, we generated shocks by evolving 
general relativistic voids with = 50, eout('^i)/c^ = 10^, uj = 1/3, a = 10""'^°, and 
AR^ = .25 and AR^ = .5. (In addition, C = .3 and / = .01). After a shock formed, 
we calculated pb by averaging p for the first few grid points of p behind the shock. In 
addition, we found $6, e?,, and Ub similarly. Then, we plugged pb into Eqs.( |A.2^ )- 
( |A.29| ) to calculate the predicted values. The percent change is plotted in Figure 13. 
These errors are expected to decrease when the artificial viscosity is decreased, since for 
these runs, the number of grid points in the shock front is 5 — 9. 
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Table 1: Relativistic Homogeneous convergence test for Atcg < At 



AR{ti) 


(Li(g) X 106) / s{q) 


n 


R 


U 


M 


e 


20.0 


28. / — 


2.5 / — 


7.6 / — 


18. /- 


2.4 / — 


10.0 


6.9 / 2.0 


.62 / 2.0 


1.9 / 2.0 


4.4 / 2.0 


.61 / 2.0 


5.00 


1.7 / 2.0 


.15 / 2.0 


.18 / 2.0 


1.1 / 2.0 


.15 / 2.0 


2.50 


.43 / 2.0 


.038 / 2.0 


.12 / 2.0 


.29 / 1.9 


.036 / 2.1 


1.25 


.11 / 2.0 


.0095 / 2.0 


.037 / 1.7 


.084 / 1.8 


.0074 / 2.3 



Table 2: Relativistic Homogeneous convergence test for At„ < At, 



/ 


(Li(g) X 10^) / s{q) 


n 


R 


U 


M 


e 


.02 


21. / 


2.2 / 


6.7 / 


16. / 


2.2 / - 


.01 


5.1 / 2.2 


.46 / 2.3 


1.4 / 2.3 


3.3 / 2.3 


.45 / 2.3 


.005 


1.1 / 2.2 


.10 / 2.2 


.34 / 2.0 


.79 / 2.1 


.094 / 2.3 


.0025 


.27 / 2.0 


.024 / 2.1 


.10 / 1.8 


.23 / 1.8 


.016 / 2.6 



Table 3: Collapse Times for a Superhorizon-Sized relativistic void 25 times the Hubble 
radius 







Atpeak area 


a = W 


a = 10-'^ 


a = IQ-''-" 


6 


1.88 


.33* 


2.05 


.047 


2.22 


.0076 


.056 


18 


1.79 


.34* 


1.99 


.047 


2.17 


.0076 


.041 


50 


1.73 


.40* 


1.94 


.049* 


2.13 


.0076 


.020 


100 




1.91 


.056* 


2.11 


.0077 


.010 


200 




1.88 


.067* 


2.08 


.0079* 


.0052 


$in 


A^colIapseZ-^^out (^i) 






.44 


.079 


.014 




.6q;~-^ 


.33 


.042 


.0052 
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Table 4: Collapse Times versus AR 





R-w ^inncr 


ar/AR^ 


A^collapsc/-^out (^i) 




1 


13 


2.05 


.047 


.29 


.5 


7 


2.13 


.044 


.30 


.25 


3 


2.20 


.034 


.31 



Table 5: Shock Location versus A7?„ for an uncompensated SR void 



AR^, 


Rsh 


Pb 




.01 


18.4 


.013 


3.59 


.02 


18.4 


.014 


3.60 


.05 


18.5 


.013 


3.57 


.1 


18.5 


.011 


3.48 


.2 


18.6 


.0085 


3.33 


Shock Tube 


18.4 


.012 


3.62 
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